function gauge_n2o_intercomparison_lifetime_get_geoschem

  fnames=file_search('/data/shared/GAUGE/N2O/GEOS-Chem/netcdf_lifetime_N2O/*.nc')

  
  time=!null
  life=!null
  q=!null
  foreach file, fnames do begin

    time=[time, julday(indgen(12) + 1, 15, mr_ncget('date', filename=file))]
    loss=mr_ncget('loss', filename=file)
    burden=mr_ncget('burden', filename=file)
    loss=[loss[0], loss[1:11] - loss[0:10]]
    burden=[burden[0], burden[1:11] - burden[0:10]]
    life=[life, burden/(loss*24.*3.*365.25)]

    qGrid=reform((mr_ncget('Emissions', filename=file))[0, *, *, *])
    area=mr_ncget('sf_area', filename=file)
    qTemp=fltarr(n_elements(qGrid[0, 0, *]))
    for ti=0, n_elements(qGrid[0, 0, *])-1 do begin
      qTemp[ti]=total(qGrid[*, *, ti]*area)/!mr_const.avo*!mr_const.mn2o/1.e12*365.25*24.*3600.
    endfor
    q=[q, qTemp]
    
  endforeach
  
  return, {life:life, time:time, q:q}
  
end


pro gauge_n2o_intercomparison_lifetime

  tomcattime=mr_mozdate(mr_ncget('date', $
    filename='/data/shared/GAUGE/N2O/TOMCAT/LIFETIME_CH4_N2O_TOMCAT3D_yrs_global_1985_2011.nc'))
  tomcatLife=mr_ncget('N2O_LIFETIME', $
    filename='/data/shared/GAUGE/N2O/TOMCAT/LIFETIME_CH4_N2O_TOMCAT3D_yrs_global_1985_2011.nc')
  tomcatLon=mr_ncget('lon', $
    filename='/data/shared/GAUGE/N2O/TOMCAT/FLUX_N2O_TOMCAT3D_ppb_global_1985_2011.nc')
  tomcatLat=mr_ncget('lat', $
    filename='/data/shared/GAUGE/N2O/TOMCAT/FLUX_N2O_TOMCAT3D_ppb_global_1985_2011.nc')
  tomcatQgrid=mr_ncget('flux_n2o', $
    filename='/data/shared/GAUGE/N2O/TOMCAT/FLUX_N2O_TOMCAT3D_ppb_global_1985_2011.nc')
  area=mr_areagrid(tomcatlon, tomcatlat)
  tomcatQ=fltarr(n_elements(tomcatTime))
  for ti=0, n_elements(tomcatTime)-1 do begin
    tomcatQ[ti]=total(tomcatQgrid[*, *, ti]*area*365.25*24.*3600.)/1.e9
  endfor

  geosOut=gauge_N2O_intercomparison_lifetime_get_geoschem()
  geosTime=geosOut.time
  geosLife=geosOut.life
  geosQ=geosOut.Q
  
  mzOut=gauge_mozart_read_history('N2O', 'N2O', 'N2O_LIFETIME_avrg', 2004)
  mzLife=mzOut.var[0, 0, *]/365.25/24./3600.
  mzTime=mr_mozdate(mzOut.date)
  mzOut=gauge_mozart_read_history('N2O', 'N2O', 'N2O_SRF_EMIS_avrg', 2004)
  mzQgrid=mzOut.var
  area=mr_areagrid(mzOut.lon, mzOut.lat)
  mzQ=fltarr(n_elements(mzTime))
  for ti=0, n_elements(mzTime)-1 do begin
    mzQ[ti]=total(mzQgrid[*, *, ti]*area*365.25*24.*3600.)/1.e9
  endfor
  
  pl=hash()
  pl[0]=plot(mzTime, mzQ, color='blue', $
    name=strcompress('MOZART, ' + string(mean(mzQ), format='(I)') + ' Tg yr$^{-1}$'), xtickunits='YEARS')
  pl[1]=plot(geosTime, geosQ, color='red', $
    name=strcompress('GEOS-Chem, ' + string(mean(geosQ), format='(I)') + ' Tg yr$^{-1}$'), /over, /current)
  pl[2]=plot(tomcatTime, tomcatQ, color='purple', $
    name=strcompress('TOMCAT, ' + string(mean(TOMCATQ), format='(I)') + ' Tg yr$^{-1}$'), /over, /current)
  pl[0].xrange=[julday(1, 1, 1996), julday(1, 1, 2013)]
  leg=legend(/auto_text_color, target=pl, shadow=0, position=[0.5, 0.95])
  pl[0].ytitle='N$_2$O global flux (Tg/yr)'
  pl[0].save, strcompress('/data/shared/GAUGE/N2O/Intercomparison/N2O_flux.png', /remove_all), $
    resolution=100
  pl[0].close

  pl=hash()
  pl[0]=plot(mzTime, mzLife, color='blue', $
    name=strcompress('MOZART, ' + string(mean(mzLife), format='(I)') + ' yr'), xtickunits='YEARS')
  pl[1]=plot(geosTime, geosLife, color='red', $
    name=strcompress('GEOS-Chem, ' + string(mean(geosLife), format='(I)') + ' yr'), /over, /current)
  pl[2]=plot(tomcatTime, tomcatLife, color='purple', $
    name=strcompress('TOMCAT, ' + string(mean(tomcatLife), format='(I)') + ' yr'), /over, /current)
  pl[0].xrange=[julday(1, 1, 1996), julday(1, 1, 2013)]
  leg=legend(/auto_text_color, target=pl, shadow=0, position=[0.5, 0.35])
  pl[0].ytitle='N$_2$O lifetime (years)'
  pl[0].save, strcompress('/data/shared/GAUGE/N2O/Intercomparison/N2O_lifetime.png', /remove_all), $
    resolution=100
  pl[0].close


end